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1.  Introduction.  We  examine  a  scalar  input-output  system  that  models  a  bound¬ 
ary  feedback  scheme  for  the  damping  of  torsional  vibrations  in  a  (not  necessarily  uni¬ 
form)  cylindrical  rod  of  circular  cross  section,  consisting  of  a  linear  viscoelastic  material. 
The  open  loop  transfer  function  for  the  system  is  irrational,  and  we  study  the  problem  of 
approximating  some  ideal  compensator  by  a  proper  rational  one.  The  approach  is  in  the 
spirit  of  [12],  where  bending  vibrations  in  an  Euler-Bernoulli  beam  with  Kelvin- Voigt 
damping  were  studied  (and  where  a  discussion  of  potential  applications  of  the  method 
is  found).  In  particular,  a  compensator  derived  from  the  full  distributed  parameter 
model  is  approximated.  The  general  method  is  designed  to  deal  with  a  wide  range  of 
viscoelastic  materials  and  structures.  Here  we  treat  the  full  range  of  linear  visco‘el astir 
constitutive  relations,  but  we  examine  only  one  structure,  namely,  the  case  of  torsional 
vibrations  (i.e.  the  viscoelastic  wave  equation)  with  actuator  and  seinsor  collocated  ar 
one  end  of  the  rod.  An  analogous  study  for  bending  motion  in  beams  will  follow  in  a 
separate  paper  [10] ;  for  non-collocated  sensor  and  actuator,  new  issues  arise  that  will 
not  be  addressed  here,  see  [12].  We  remark  that  in  addition  to  the  work  in  [12],  the 
problem  of  robust  controller  design  for  an  Euler-Bernoulli  beam  with  strong  damping, 
namely  Kelvin- Voigt  damping,  has  been  considered  in  several  papers;  in  particular,  u  c 
cite  [11],  [3]  and  [4]. 

In  the  particular  (collocated)  cases  that  we  discuss,  the  open-loop  transfer  fimction.s 
P{s)  have  no  zeros  or  poles  in  the  open  right  half-plane,  and  no  zeros  or  poles  on  th(' 
imaginary  axis  apart  from  a  pole  at  zero  and  a  (fractional  order)  zero  at  infinity.  The 
transfer  function  of  a  typical  “optimal”  compensator  will  be  some  rational  function 
divided  by  P{s);  since  the  latter  function  is  irrational,  the  compensator  obtained  in 
this  way  will  be  irrational  too.  In  order  to  get  a  compensator  that  can  be  physically 
implemented,  one  has  to  approximate  P{s)  by  a  rational  function,  preferably  of  low 
degree. 

As  a  model  problem  illustrating  the  type  of  approximation  we  will  make,  consider 
the  equation 

Wit{x,t) -\-2aWi(x^t)  +  a^w{xyt)  =  Wxx{3:,t) 

(0  <  .T  <  1,0  <  f  <  00,G  >  0) 

with  boundary  and  initial  conditions 

^Ga:(0,t)  =  0,u;i(l,t)  =  uit)  =  input,  u»(n;,0)  =  iut{x^0)  =  0 

and  output  This  is  a  wave  equation  with  dispersion  and  viscous  damping,  and 

it  is  not  formally  covered  by  our  results,  but  it  is  quite  similar  in  nature.  Solution  of 
the  system  for  the  Laplace  transform 

fOO  _ 

u)(a?, s)=  /  e  ^^w{x^t)dt 
Jo 

leads  to  the  formula  ui(l,  5)  =  P(s)u(s),  where  the  transfer  function  is 

m  P(.)  -  cosh(5  +  g)  ^  1  fr  1  +  (s  +  a)Vg| 

(5  -i-  a)  sinh(s  -Ha)  (s  -H  a)^  1  -f  (5  -f  a)^/?^^ 
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with  ^  =  (2fc  —  l)?r/2,  Tjk  =  kw.  We  attempt  to  approximate  P  with  a  finite  product 


(3) 


Pn{s)  s 


1 

(s  +  a)2 


n 


it=:l 


l  +  {s  +  a)y4- 


As  we  shall  see  in  Section  2,  this  approximation  should  be  done  in  such  a  way  that, 
e.g.,  {P{s)/Pi^{s)  -  l)f{s)  tends  to  zero  uniformly  in  the  right  half-plane;  the  function 
T  is  a  certain  weight  function  (a  design  parameter)  that  has  a  zero  at  infinity. 

In  a  uniform  or  nonimiform  viscoelastic  rod  model  one  gets  a  similar  transfer  func¬ 
tion.  The  main  difference  is  that  the  factor  (s  +  a)  is  replaced  by  a  (possibly)  tran.sccii- 
dental  function  ^(s)  which  reflects  the  viscoelastic  memory  kernel  for  the  material. 

The  geometry  of  the  problem  remains  the  same  for  all  different  functions  /i(s),  so 
the  same  partial  product  scheme  is  used  for  all  materials.  We  prove  convergence  of  these 
partial  products  in  regions  of  the  plane  that  include  the  range  of  j3{s).  The  dotail(?d 
mechanical  properties  of  the  viscoelastic  material  enter  solely  through  0(s)  and  will 
influence  the  number  of  factors  needed  for  the  desired  degree  of  accuracy  in  the  partial 
product. 

The  product  expansion  that  we  use  is,  in  principle,  well-known.  The  functions 
that  need  to  be  approximated  have  alternating  zeros  and  poles  on  straight  lines  in  the 
complex  plane.  In  particular,  in  the  case  of  the  rod  these  zeros  and  poles  lie  on  the 
imaginary  axis.  According  to  the  Herrnite-Biehler  theorem,  functions  of  this  type  can 
be  represented  as  infinite  products  of  linear  fractional  transformations.  However,  the 
Herrnite-Biehler  theorem  itself  does  not  give  a  bound  on  the  behavior  of  the  error  terms 
at  infinity,  and  a  substantial  part  of  this  work  is  devoted  to  the  development  of  a  precise 
description  of  the  asymptotic  error  at  infinity  (Section  4).  In  particular,  we  show  that 
one  gets  a  better  approximation  at  infinity  if  instead  of  simply  truncating  the  infinite 
product  expansion,  we  also  multiply  the  truncated  product  by  an  additional  irrational 
factor  (sec  the  discussion  preceding  equation  (31)  in  Section  4). 

To  complete  the  construction  of  a  rational  approximation  of  P  for  models  wher(' 
is  irrational  we  have  to  approximate  linear  fractional  transformations  of  .d^{s)  hy 
linear  fractional  transformations  of  $.  In  Section  5  we  suggest  a  very  simple,  low  order 
approximation  which  in  many  cases  leads  to  satisfactory  results.  The  idea  behind  thi.s 
approximation  is  to  separate  the  dynamic  modes  from  the  creep  modes,  and  to  ignore 
the  latter  type  of  modes.  This  seems  to  work  quite  well  when  the  internal  structural 
damping  is’  small. 

Finally,  in  Section  6  we  discuss  the  extra  irrational  convergence  factor  that  was 
introduced  in  Section  4,  formula  (31).  In  particular,  we  show  with  as^anptotic  estimatc-s 
and  numerical  experiments  that  we  get  reasonably  good  results  by  simply  ignoring  the 
extra  factor  completely,  provided  that  the  design  parameters  S  and  T  are  propcily 
chosen. 

The  authors  thank  Professor  Joseph  A.  Ball  for  several  helpful  discussions. 
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2.  Weighted  Rational  Approximations  of  the  Optimal  Compensator.  As 
we  shall  see  in  Section  3,  the  transfer  functions  that  we  deal  with  here  can  be  expressed 
as  the  quotient  of  two  outer  iy®®-functions.  More  precisely,  for  the  type  of  rod  problem 
treated  here,  they  are  continuous  and  have  no  zeros  or  poles  in  the  closed  riglit  half¬ 
plane,  except  for  a  zero  at  infinity  (not  of  integer  order,  in  general)  and  a  pole  of 
order  two  at  zero.  (In  the  essentially  identically  situation  where  the  left  end  of  tlic 
rod  is  fixed  rather  than  free,  the  transfer  function  is  regular  at  zero.)  In  particular, 
this  means  that  by  the  appropriate  choice  of  compensator,  one  can  achieve  almost  any 
sensitivity  function  5  =  (1  -f-  The  only  constraints  that  have  to  be  imposed  on 

S  are  some  growth  conditions  at  infinity,  and  possibly  at  zero. 

Within  this  class  of  feasible  sensitivity  functions  5,  one  would  like  to  find  an  optimal 
one.  There  are  several  different  methods  available,  including  mixed  sensitivity 
optimization  and  P^-optimization  with  constraints.  It  may  happen,  as  in  [18], 
that  S  itself  is  already  a  suboptimal  approximation  of  an  ideal  sensitivity  function 
that  fails  to  satisfy  the  stability  conditions  (4)  below.  Here  we  shall  not  go  into  this 
optimization  procedure.  Instead  we  assume  that  one  has  found,  in  one  way  or  another, 
a  sensitivity  function  S  that  one  would  like  the  closed  loop  system  to  have. 

As  we  mentioned  above,  not  every  S  is  permissible,  because  the  open  loop  transfer 
function  P  has  a  zero  at  infinity,  and  possibly  a  pole  at  zero.  This  can  be  seen  in  many 
ways,  one  of  which  is  the  following.  The  closed-loop  system  is  alwa3^s  required  to  be 
stable,  i.e.,  the  four  closed-loop  transfer  functions 


S  = 

f  = 

CS^ 

= 


1 

14- PC’ 

1  -  5  =  PCS  = 
C 

14- PC’ 

P 

14- PC’ 


PC 

14- PC’ 


should  all  belong  to  P®®,  and  their  P®®-norms  should  not  exceed  certain  numbers  that 
can  be  given  a  priori.  (These  numbers  depend  on  the  sizes  of  the  sensor  and  actuator 
noises,  etc.)  In  particular,  we  have 


(4) 


IICSII//00  <  Kc.  ||Ml|/,~  <  A>. 


for  some  constants  Kc  sii^d  Kp.  For  our  scalar  system  the  former  of  these  conditions  is 
equivaleut  to 

(5)  |1- 5(5)1  =  If  (^)|<Ac|P(s)|, 
and  the  latter  is  equivalent  to 

(6)  |5(5)|  <  A>|P(5)r*. 
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Note  that  (5)  gives 

(7)  5(5)  =  1  +  0(|P(5)|)  as  5  00, 

which  combined  with  (4)  implies 

limsupjC(s)|  <  Kc, 

that  is,  C  has  to  be  proper.  The  restriction  (6)  comes  into  play  mainly  in  those  cas{'s 
where  P  has  a  pole  at  zero,  since  it  implies  that 

(8)  S{s)  =  0(|P(s)|-i)  as  s  0. 

The  conditions  (4)  imply  stability.  To  see  this,  observe  that 

Sf  =  {CS){PS)] 

hence,  for  each  s, 

|5(^)l(|5(s)|  -  1)  <  KcKp, 


and 

|^(s)l  <  1/2  +  \/l/4  +  Kcl^p. 

A  similar  estimate  holds  for  |T(s)|. 

In  the  sequel  we  assume  that  the  “optimal”  sensitivity  function  S  that  we  would 
like  to  implement  is  feasible  in  the  sense  that  it  satisfies  (4). 

Once  a  sensitivity  function  S  is  given,  one  may  compute  the  corresponding  optimal 
compensator  from  the  formula 

f 

SP  SP' 

There  is  only  one  thing  wrong  with  this  C;  in  general  it  will  not  be  rational,  due  to  the 
fact  that  P  is  not  rational,  and  it  cannot  be  implemented  exactly  by  means  of  standard 
circuits.  This  leads  to  the  main  point  of  this  paper.  Suppose  that  we  have  a  given  S 
satisfying  (4).  How  can  we  find  a  rational  compensator  such  that  the  sensitivity 
function  Sj^  corresponding  to  this  compensator  will  be  reasonably  close  to  the  ideal  one 
5,  and  will  satisfy  the  analogue  of  (4),  namely, 

(9)  <  Kc,  <  -fvp, 

as  well?  We  shall  assume  that  the  optimal  S  (hence  the  optimal  f)  is  rational.  Thus, 
it  is  really  a  question  of  approximating  P  by  a  rational  function  Pv,  and  setting  Ca^  = 

f/iSP^). 
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There  are  several  ways  that  one  may  interpret  the  words  “reasonably  close”  in 
the  preceding  paragraph,  and  to  each  different  way  corresponds  a  different  measure  of 
how  good  the  approximation  is.  Perhaps  one  of  the  most  natural  measurements  of  tlio 
goodness  of  an  approximating  compensator  is  the  smallness  of  either  the  i^^-norni 
or  the  H^-norm  of  -  §.  These  norms  can  be  written  in  the  form 


\\Sr^  -  S||//p 


T  +  CatP  l  +  CP\\jj, 

(  1)  ^  II 


where  p  =  2  or  p  =  oo,  and 

r^  =  C^/C  =  P/P^. 

If  we  ignore  higher  order  terms,  then  this  becomes  \\{Tf^r  —  l)ST||y/7..  Thus,  in  this  case 
the  problem  could  be  interpreted  as  finding,  for  each  N,  a  controller  C^i  of  order  A^ 
say,  that  minimizes 

In  other  words,  we  get  a  weighted  approximation  problem,  where  the^weight  is  ST. 
Because  of  (7),  this  weight  function  will  have  a  zero  at  infinity,  and,  if  P  has  a  pole  at 
zero,  then,  because  of  (8),  it  will  also  have  a  zero  at  zero. 

In  the  discussion  above  we  have  ignored  all  higher  order  terms  and  the  stalnlity 
bounds  (9).  When  these  are  taken  into  account  one  gets  some  additional  constraints. 
The  exact  expression  for  5^  —  5  is 

Spf  —  S  =  {I  —  Tyv)T5Ar, 


where  Sn  can  be  written  in  the  form 


(10) 


Sn  = 


S 

1^{TN-I)t' 


Thus,  we  have  the  more  precise  estimate  (for  p  =  2  or  p  =  oo) 

||S«  -  5||;„  <  |I(r^  -  l)5f  ||(l  +  (r;.  -  . 

To  make  this  norm  small  it  suffices  to  minimize  the  P^-norm  of  {tn  —  l)ST  suljjcct  to 
the  constraint  that  for  some  constant  a  <  1, 

(11)  -  I)r||^„  <  a  <  1. 

This  is  the  same  problem  as  before,  except  for  the  additional  constraint  (11).  Of  course, 
one  could  also  use  the  fact  that 


.  ||{-.-i)^?lL,.<il(-^-i)?LI|5|!„„. 
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and  minimize  the  H^-novm  of  (tat  —  1)T  instead  of  the  iif^’-nonn  of  (ta*  —  1)5?  since' 
the  ^“-norm  of  5  is  finite. 

Let  us  still  take  a  closer  look  at  the  stability  bounds  (9).  Because  of  (4),  (10)  and 
(11),  we  have 


||(l  +  (tjv  -  l)f)  '1^^^ 


and 


<—■ 

//«=  1  —  a 


Thus,  our  problem  becomes  the  one  of  minimizing  the  of  either  [r.\-  -  1)57' 

oi'  ('Tftr  —  1)T  subject  to  the  constraint  (11),  combined  with  cither 


(12) 

VI 

i: 

or 

(13) 

where  /?  >  1  is  some  constant.  Note  that  (12)  can  be  much  more  conservative  than  ( 13) 
due  to  the  fact  that  whenever  the  function  CS  has  a  zero  at  infinity  it  i.s  i>ossil)lc  lo 
allow  Tat  to  be  quite  large  or  even  infinite  at  infinity;  this  is  true  whenever  C  Is  strict  Im¬ 
proper. 

If  one  chooses  some  other  criterion  for  goodness  of  the  approximation,  then  oin: 
gets  some  other  weighted  approximation  problem,  either  in  or  in  with  addi¬ 
tional  constraints.  See,  e.g.,  [1].  These  weights  will  always  have  a  zero  at  infiuitm', 
and  sometimes  a  zero  at  zero,  and  ^^“-constraints  of  the  type  (11)  and  (13)  are  al¬ 
ways  present,  either  explicitly  or  implicit,  since  they  are  intimately  connected  to  the 
stability  of  the  feedback  system. 

The  method  that  we  propose  will  not  be  optimal  with  respect  to  any  of  thesf' 
minimization  problems;  indeed,  the  approximate  plant  Pjv  will  not  be  charactcrizcHl  as 
a  minimizer.  The  method  will,  however,  give  fairly  good  results  with  respect  to  all 
measures  of  closeness  of  the  tj^pe  described.  Specificalty,  we  shall  construct  rational 
approximates  Fpf  (of  orders  that  will  increase  linearly  in  N)  in  such  a  way  that 

•  ri\f  =  Ci^lC  =  P/Pn  1  uniformly  on  compact  subsets  of  the  right  half-plane, 
and 

•  limsupjv^^||rArC5||i?~  <oo. 

It  is  easy  to  show  that,  for  the  plants  that  we  consider,  and  for  any  compensator 
C  satisfying  the  stability  bounds  (4),  we  will  then  have  \\{Ti^r  —  1)5T1|j/p  ^  0  and 
IK^Tv  -  0  for  p  =  2  and  p  =  oo  (note  that  C{s)S($)  C{5)  and  T(.s)  ~ 

C{s)P{s)  as  l^l  — ^  oo).  In  particular,  (11)  will  be  satisfied  for  large  values  of  N  with 
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arbitrarily  small  a,  and  (13)  will  hold.  Of  course,  our  compensators  will  not  have  the 
lowest  possible  order  corresponding  to  a  given  accuracy,  but  the  order  may  be  further 
reduced  by  means  of  some  standard  order  reduction  scheme. 

As  we  shall  see,  we  will  not  quite  achieve  the  goal  set  forth.  First,  if  P^ia)  i.s  not 
rational,  our  method  of  constructing  Pn  (Section  5)  leaves  a  residual  error  due  to  some 
continuous  creep  modes  (singularities  other  than  poles  on  the  negative  real  iixis)  that 
we  ignore.  Although  this  residual  error  does  not  tend  to  zero  as  TV  oo,  it  will  hv 
quite  small  in  the  case  where  the  material  is  “nearly”  elastic. 

Second,  our  best  approximations  involve  multiplication  of  the  basic  partial  product 
by  a  new  factor,  irrational  in  and  unbounded  as  s  oo.  In  Section  6  we  analyze  the 
effect  of  dropping  or  simplifying  this  factor;  the  same  type  of  convergence  (A^  oc) 
holds,  but  the  convergence  rates  are  weakened.  A  detailed  investigation  of  low  orrl(}r 
rational  approximates  for  this  factor  remains  to  be  done. 

At  each  stage  of  the  approximation,  we  introdiice  a  new  relative  error,  in  Section 
4,  in  Section  5  and  in  Section  6,  and  each  of  these  sections  will  close  with  some 
plots  of  |r^^(ia;)  -  1|  for  some  model  kernels.  The  cumulative  relative  error  resulting 
from,  making  all  the  approximations  is  then 

We  now  proceed  to  give  a  detailed  description  of  our  method. 
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3.  The  Physical  Models  and  their  Transfer  Functions.  We  stud}-  toisioiial 
vibrations  in  a  cylindrical  rod  of  circular  cross  section  consisting  of  a  viscoelastic  niafe- 
rial^  for  which  one  modifies  the  usual  equations  of  motion  from  linear  elasticity  theory 
by  permitting  the  stress  at  any  point  to  depend  on  the  history  of  the  strain  rate  at  that 
point.  The  memory  dependence  is  governed  by  a  relaxation  function  A{t)  =  E  o(f). 
where  E  >  0  and  a{t)  is  completely  monotonic  with  0  =  a(oo)  <  ft(04-)  <  oo  and 

/•oo 

/  dt  <  oo  for  some  6  >  0. 

Jo 

As  an  extreme  case  we  also  include  a{t)  =  a6{t)  where  ft  is  a  positive  constant  and 
<5(t)  is  the  delta  function  (Kelvin- Voigt  damping).  Note  that  the  case  of  pure  elasticity, 
a{t)  =  0,  is  excluded. 

Reasoning  as  in  [9],  we  arrive  at  the  equation 

(14)  p(x)Wi{x^t)  =  [  A(t  —  r)Wj;j;{x^r)dr  (0  <  .r  <  1,  t  >  0). 

Jo 

Here  (after  scaling)  w{x^t)  represents  the  torsion  at  position  x  along  the  rod  at  riiiu'  ^ 
The  density  function  p{x)  is  assumed  to  be  smooth  (C^)  and  strictly  positi^'o  in  [0. 1], 
We  assume  that  the  rod  is  at  rest  for  t  <  0  and  that  its  left  end  is  stress-free.  At  the 
right  end  we  give  a  stress  input  u{t)  via  an  applied  torque  and  measure  as  nii 

output.  Thus, 


ty(a;,0)  =0  (0  <  a:  <  1), 

(15)  <7(0, t)  =  r)w^{0,  r)<ir  =  0  {t  >  0), 

(T(l,f)  =  -^  /  A(t  —  r)wi(l,r)cl7'  =  ti{t)  {t  >  0). 
dt  J 0 

We  shall  be  working  in  the  frequency  domain;  in  other  words,  we  study  the  La]:>lac(' 
transforms  of  the  quantities  of  interest,  denoted  by  0(5),  u;(*t,5),  etc.  We  introduce 


a{s)  =  {sA{s)y^‘^  =  (F  +  sa{s)Y^^  {seC\  (-00, -^ol), 


P{s)  =  s/a{$). 

Here  and  below  is  defined  so  that  >  0  when  2  >  0;  <5o  is  the  largest  number 
such  that  60  <  <5  and  sA(s)  >  0  on  (-5o,0)-  After  applying  the  Laplace  transform  tc^ 
equations  (14)-(15),  we  get 

a^{s)wxx{^,  s)  $^p{x)w{x,  s)  =  0  (0  <  a;  <  1) 

(16)  a^(s)u;a.(0,  s)  =  0  *  (free  end  condition) 

o:^(s)u;a;(l,  5)  =  u(s)  =  boundary  stress  input. 

When  p{x)  =  1  the  differential  equation  and  the  first  boundary  condition  in  (1C) 
have  the  general  solution 

w{x^s)  =  0(5)  cosh )0a:  with  (J  =  P{s)^ 
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so  we  get  that  ^(l,s)  =  P{s)u{s),  where  the  transfer  function  P  is  given 
(17)  P(s)  =  =  ^coth/3{s). 


Notice  that  /(/^)  has  the  product  expansion 


/{/?)  =  coth/7=in 


l  +  (/3/6)^ 
l  +  WVkf’ 


with  zeros  and  poles,  respectively,  at  the  points 


P  =  ±ffci  =  and  /?  =  r/o  =  0,  P  =  ±%i  =  ±fciri  (fc  =  1, 2, . .  .)■ 

For  the  case  of  general  densities  below,  we  have 

Proposition  3.1.  The  transfer  function  for  the  rod  problem  has  the  foim 


(18) 


where  /(/?)  has  the  factorization 


Po 


Here  0  =  tjo  <  Ci  <  ’?i  <&<%••  • 

(19)  +  0(l/fc)  and  nk  =  —  +  O(lA-)  {k  -  oo), 

2pi  pi 

and  pi  =  Jq  p{i'y^^dx. 

Proof  Fix  s  (and  hence  /?  =  /?{s))  and  let  yix,z)  be  the  solution  of 

(20)  2/"  +  z^p{x)y  =  0,  2/(0)  =  1,  2/(0)  =  0. 


As  in  the  uniform  case,  we  can  then  compute 


w{l^s) 


2/(1, 0) 

52  2/(1,  i/?) 


u{s) 


=  Pis)u{s). 


To  develop  the  product  expansion  claimed  for  /(/?)  =  /?2/(l.  >-P)<  "'c  ajjply  [13, 

Tlim.  1,  p.  308,  and  Thin.  4,  p.  316]  to  the  function 


g{z)  =  zy{l,z)ly\l,z). 


We  establish  that 
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(i)  y{l^z)  and  y\l,z)lz^  are  entire  and 

(21)  ~~  Ja 


(ii)  y{\^z)  and  y’{\^z)lz  have  no  common  zeros, 

(iii)  if  m  and  n  are  real  and  z  is  nonreal,  then  my{l,z)  +  ny\l,  z)lz  7^  0, 
and 

(iv)  g{z)  is  increasing  on  the  real  axis  (on  the  intervals  between  its  poles). 
Once  (i)-{iv)  are  proved,  the  results  from  [13]  yield  the  expansion 


g{2)  =  c 


z-ap 

z 


n 

fcGZ\{0} 


(1  -  z/ak) 
(1  -  z/h)  ’ 


where  c>  0  and  (with  60  =  0),  h  <  <  bk^i  for  eveiy  integer  k.  Since  y(lA))  ~  1. 

relation  (21)  then  shows  that  oqC  =1/  Sq  p{^)  da;. 

From  (20)  it  is  clear  that,  aside  from  bp,  the  a^s  and  b^s  occur  in  posit ive-negati\e 
pairs  so  the  desired  product  expansion  for  /(/?)  follows.  Finally,  the  asyinj)- 

totic  distribution  of  the  and  t}j,  follows  from  Sturm-Liouville  theory  [2,  p.  328).  [5,  p. 
414]. 

Assertions  (i)  and  (ii)  follow  from  the  elementary  theory  of  differential  equations. 
In  particular,  for  the  an^yticity  of  y'{liZ)fz^  at  0  and  (21),  differentiate  (20)  with 
respect  to  z  to  get 

z)  +  z^p{x)y.Xx,  z)  =  -2zy{x,  z),  j/,(0,  z)  =  2/',(0,  i )  =  0. 


(We  denote  dy/dx  by  2/^,  dyjdz  by  t/.,  etc.;  all  the  derivatives  that  will  appe^nr  nu' 
continuous  in  {x,z),)  Letting  ^  0,  we  get  that  yt{x,0)  =  0,  0  <  ;i‘  <  1.  so 

lim,_,o2/(l»'2^)/^  =  2/Ul»0)  =  0.  Differentiating  again  in  2;,  we  get 

z)  =  -4zp{x)y^{x,  z)  -  2p(x)y[x,  z). 


2/«(0,^)  =  s4,(0,z)  =  0. 

Since  ?/(a;,0)  =  1,  we  can  let  2:  0  to  obtain  24^(0;,  0)  =  — 2p(a;),  so 

lini7/'(l,z)/z“  =  p(x)d.r, 


as  asserted. 

Assertion  (iv)  follows  from  Sturm-Liouville  theory.  From  [2,  p.  312,  (22)],  witli 
fl(z)  =  arctan(2/(l,z)/y'(l,z)), 

^  =  /9(1)  sin^  9  +  cos^  5  >  0,  z  €  R, 
az 

so  g{z)  is  a  product  of  two  increasing  functions  between  the  zeros  of  ?/^{l,  z)  on  R. 
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Finally,  for  assertion  (iii),  suppose  there  is  a  ^  G  C\R  such  that  zy{l,z)ly'(\,z)  = 
r  G  R\  {0}.  (By  standard  boundary  value  theory  y{l^z)  =  0  and  =  0  arc 

impossible.)  Multiply  (20)  by  y  and  integrate  by  parts  to  get 

(22)  ”  ”|2/(^)I"* 

Now  set  z  =  p  +  iu  and  separate  (22)  into  real  and  imaginary  parts: 

0  =  I/^)  p{x}\y{x}\‘^  dx  -  |j/(x)p  dx  -  ^13/(1)!" 

0  =  2pi'  f  p(x)|i/(x)P  dx  -  ^|j/(l)P- 
Jo  T 

Since  z/  9^  0,  we  can  cancel  u  in  the  second  equation  and  substitute  for  |2/(l)r“/r  in  the 
first  to  obtain 

0  =  f  p(^')b(-'*^)P 

Jo  JQ 

which  is  impossible.  D 

We  remark  that  the  asymptotic  error  0{k^^)  in  formula  (19)  cannot  be  improved 
in  general  [15,  Section  4.3]. 

Next,  recall  from  the  Introduction  that  we  propose  to  approximate  P(s)  ~ 
6‘"^/?(5)/(/5(^))  by  rational  functions  in  two  steps,  where  in  the  first  step  wc  approxi¬ 
mate  Pf{l3)  by  a  rational  function  of  p.  The  domain  of  s  that  is  of  concern  is  the  right 
half-plane  >  0.  Thus,  we  need  a  rational  approximation  of  /  that  is  valid  in  tlie 
image  of  the  right  half-plane  under  /?,  i.e,,  the  domain  of  our  approximation  is 

n  =  {/?(s)15?s>0}. 


We  complete  this  Section  by  describing  this  set. 

The  set  11  depends  heavily  on  the  size  of  the  function  .4  and  its  derivative  in  the 
neighborhood  of  zero. 

Lemma  3.2.  Under  the  hypotheses  of  this  section, 

L  ra;/?(0)  =  0, 

(h)  ^p{s)  for  Us>0,s^0,  . 

(cj  p{s)  =  Pis), 

(d)  ^p{iuj)  -4  00  a5  w  00, 

(e)  ^p{iJ)  >  ^p(iuj)  for  u>0, 

(f)  {arg(s)<5r/4}cn  =  {/5(s)l»«>0,  s  5^  0}  C  {  arg{s)  <  7r/2  }, 

(g)  liminf5_eo|«I”‘''^l/5(^)l  >  0  uniformly  in  Ss  >  0, 

(k)  P{s)  ~  if  >1(0+)  <  00  and  I3(s)  =  o{s)  if  .4(0+)  =  00  as 

5  — >  00  uniformly  in^s  >  0 . 

2.  IfA{0+)  <  oo,  then  ^p{iuj)  =  oi^Piiu;))  cs  w  oo.  Thus,  11  is  not  contained 
in  any  sector  of  the  form  {  arg(s)  <  7r/2  -  e}  for  any  e  >  0.  In  addition, 
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(23) 


^P(iij) 


a;  00, 


and  satisfies  the  regularity  condition 


(24) 


lim  inf 


^P{iru) 


w-ool<r  3^^(icx;) 


>  for  each  c  G  (0, 1), 

>  1. 


3.  If  A'{0-h)  >  -oo  (hence  ^4(0+)  <  oo),  then  ^/?(iu,’)  -  tt;-4(0+)"^/-  0 

and  ^p{iu)  '-“>l'(0+)A(0+)“^/^  as  u  oo.  7/ .4'(0+)  =  -cx:.  then 
^P{iu))  — +  00  and  UP{iuj)  — >  oo  as  u  — >  oo.  hi  particular,  in  all  cases. 

lim  infu,_»oo  >  0. 

Proof  1.  We  have 


(25) 


E  +  5ft(5)  ’ 


where  [9,  (2.14)]  m{s)  >  0  (9F^5  >  0,  5  7^  0),  S5£J2(5)  <  0  {^s  >  0,  Qf6-  7^  0),  aiul 
a  G  L^(R''‘).  Then  (la)  is  obvious,  and  a  short  calculation  shows  that  ^•s^^v5-(.s)  >  0 
(S‘5  7^  0)  so  that  |arg/?^(iw)|  <  tt,  from  which  we  see  that  p{s)  =  (s-/q-)^/-,  so  (11>)  niicl 
(Ic)  are  clear.  For  $  =  icj,  oj  >  0,  (25)  can  be  resolved  into  real  and  imaginary  parts, 
showing  that  p\itjj)  lies  in  the  second  quadrant,  so  (le)  follows.  The  fact  that  p(s)  oc 
as  s  — >  oo(3fJs  >  0)  together  with  (le)  imply  (Id).  Assertion  (If)  is  a  consequence  of 
(lb),  together  with  (la),  p{s)  — >  00  as  s  00,  (Ic),  (le)  and  the  argument  principle. 
Finally,  (Ih)  follows  from  (25)  and  [9,  (2.18)],  and  (Ig)  is  a  coascquence  of  (25)  and  tlu' 
fact  that  |a(s)|  <  a(0)  for  Rs  >  0. 

2.  We  have 


(26)  Aiiw)  = 

ice  lo; 
so 


cjjA{ioj)  =  »iA(0+)  +  0(1)  (cc;  -►  00). 


Using  this  and  (25),  we  get  that  tt  >  arg;32(ia;)  (lj  00),  and  the  first  half  of 
Part  2  follows. 

To  prove  the  regularity  condition  (24),  we  recall  that  the  completely  monotonic 
function  —A'(t)  has  a  Bernstein  representation  [17,  p.  160] 

-A\t)  =  Jf  e-^‘dv{x)  (t>0), 


w^here  i/  is  a  positive  measure,  and 


-.s.-  /•oo 


a:di^(a;) 


di/{x) 

-f 


{oj  G  R). 


From  (25)  and  (26)  we  get 


P{\u) 


fl  4,  ^'(‘‘•’)1 

A(0+)j  [  A{0A) 
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and  then  the  Taylor’s  expansion  for  (1  +  z)  yields 

Thus,  (23)  holds,  and  (24)  easily  follows  since  a;^(rr^  +  aj^)“^  is  increasing  in  u  and 
{oof  (x^  +  (ccj)2)"^  >  +  w^)"^  when  c  €  (0, 1). 

3.  This  was  proved  in  [9,  Lemmas  2.2  and  2.3].  □ 

For  particular  examples  one  can  very  precisely  describe  the  region  11.  Wc  shall 
examine  closely  the  following  model  kernels  (in  decreasing  order  of  structural  damping): 

1.  Ai(s)  =  E/s  6;  Kelvin- Voigt  damping,  where  formally  .4i  is  the  suAi  of  a 
constant  and  a  constant  times  the  unit  point  mass  at  zero. 

2.  .42(0  =  -2^  +  .42(s)  =  E/s  -|-  e(l  +  s/(5)-^  0  <  <  1,  e,  <^  > 

0,r  =  gamma  function;  a  modified  “fractional  derivative”  model  (see  [16])  of 
order  1  -  //  with  exponential  decay  as  t  — >  co. 

3.  .43(0  -  ^  Vr(M  +  1))  Sr  dr, 

.43(5)  =  E/s-\-(€6/(fis))  (1  -  (1  +  s/6y^)\  an  intermediate  model  of  order  1-//. 
with  .4(0+)  <  00  and  A'(0+)  =  -00. 

4.  A^{t)  =  F;  +  e^e"^^  .44(5)  =  jK/s  +  e/  (1  +  s/6)\  standard  linear  solid  model. 
Observe  that  the  constants  have  been  chosen  in  such  a  way  that  in  all  cases  4(.5)- 
£;/s  -  e  — ^  0  as  s  -4  0,  so  that  the  different  examples  have  the  same  low  frequency 
behavior.  The  constant  fx  is  related  to  the  behavior  of  A  near  zero,  and  6  rej^veseuts 
a  bandwidth  constant  (the  transfer  functions  do  not  differ  much  from  each  other  foi¬ 
ls]  «  <5). 

For  the  kernels  listed  above,  we  can  use  the  binomial  series  to  deduce  the  following 
precise  estimates  valid  for  large  values  of  [sj  in  SfJs  >  0.  The  proofs  are  straightforward 
and  are  left  to  the  reader. 

Lemma  3.3.  The  following  estimates  are  valid: 

1.  If  A  —  A\,  then  P(s)  =  ^  |^|  00,  SRs  >  0. 

t  //I  =  .42,  then  P{s)  =  |^.|  _  ^o,  3?.s-  >  0. 

3.  If  A  =  A3,  then  P{s)  =  (E  -¥  eS/^y^^^s  +  +  e(5///.)“V’^s^”''  + 

0(|s|^”^^)  as  jsj  ->  00,  SRs  >  0. 

4.  IfA=:A^,  thenp{s):=^{EAe6y^^h-¥le6^(E-he5y^^^AO{\s\-^)  as  |s|  oc. 

^s  >  0. 

The  estimates  in  Lemma  3.3  allow  one  to  describe  the  asymptotic  behavior  (u;  00) 
of  the  curves  /?(+iw).  In  the  following  paragraph,  C  denotes  a  positive  constant  that 
depends  on  the  model  as  well  as  the  values  of  the  parameters  f.i,E,€,S.  We  leave  it  to 
the  interested  reader  to  determine  the  value  of  C  in  each  case.  All  estimates  arc  valid 
as  tj  — >  00.  _ 

We  have  ^^(-ia;)  =  P{iuj).  For  Ai,  /?(ia;)  and  for  Ao,  P{iuj)  ^ 

Hence,  for  the  Kelvin- Voigt  model  and  for  fractional  dcl•i^’ative 
models,  FI  is  contained  in  a  proper  subsector  of  the  right  half-plane.  For  A3,  ^/3(iu;)  ~ 
C  In  particular,  for  A3, 5R/?(icj)  00  as  cj  00,  but  11  is  not  contained 

in  any  .proper  subsector  of  the  right  half-plane.  Finally,  for  the  standard  linear  solid 
model  A4,  ^p{iu})  ^  (£^  +  e(5)“^/^a>  and  9f^i0(ia;)  \e6\E  + 
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Fig.  1.  Plots  of^Q{\(jj)  (abscissa,  left)  versus  Jm^{\u)  forO  <  u  and  ofu  (abscissa,  right)  versus 

Plots  of  the  functions  P{iuj)  and  P{ioj}  for  the  different  choices  A1-.44  for  the  kc'nicl 
A  are  given  in  Figure  1.  (In  these  plots,  and  in  all  later  plots,  the  parameters  ha\-e  been 
chosen  as  follows:  =  1,  e  =  0.01,  6  =  20,  and  /i  =  0.5.)  Note  that  stronger  viscoelastic 
damping  corresponds  to  a  /?  curve  bending  more  sharply  to  the  right.  The  same  cnrv(‘ 
style  (solid  line  =  Ai,  dotted  line  =  A2y  etc.)  is  used  in  both  graphs,  as  well  as  in  later 
graphs  of  the  same  type. 
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4.  Rational  Approximation  of  the  Plant;  First  Step.  Next  we  discuss  the 
approximation  of  the  plant  in  terms  of  a  rational  function  of  and  wc  begin  witli  the 
case  of  the  rod  with  density  p{x)  =  1.  As  we  saw  in  Section  3  formula  (17),  in  this  case 
the  transfer  function  is 

P(s)  =  ^|^coth/?(s). 

If,  for  the  moment,  we  ignore  the  factor  then  we  are  left  with  the  problem  of 

getting  an  approximation  of 

/(/?(s))  =  coth/?(s). 


As  we  mentioned  earlier,  we  intend  to  do  this  approximation  in  two  steps.  First  we 
approximate  f{P{s))  in  terms  of  a  rational  function  of  /?(s),  and  then  we  approximate 
this  by  a  rational  function  of  s. 

In  the  first  step,  by  the  argument  that  we  gave  earlier,  we  need  a  rational  approxi¬ 
mation  of  the  function  /  that  is  valid  at  least  in  a  sector  of  the  type  {  ^  €  C  |  < 

7r/4},  and  the  bigger  a  sector  that  we  can  allow  the  better.  To  treat  the  case  wliere 
-.4'(0-l-)  <  00  we  actually  need  a  domain  of  approximation  11  that  is  asymptotic  to  a 
closed  half- plane  strictly  contained  in  >  0. 

There  is  an  obvious  candidate  for  a  rational  approximation.  Recall  that 


00 


f{z)  =  cothz  =  Ifz  n 

it=l 


1  +  W^ 


where 


=  nk  -  7r/2,  =  tt/c, 


and  the  convergence  is  uniform  on  compact  subsets  of  the  complex  plane,  not  containing 
any  of  the  poles  ±irjk  of  coth.  However,  we  need  convergence  in  an  unbounded  domain, 
and  we  need  fairly  explicit  estimates  on  the  error,  so  that  we  know  what  happens  after 
we  multiply  the  function  by  some  appropriate  weight  functions.  Observe  that  if  wc 
define  a  finite  approximation  fpf  by 


Mz) 


If  l±l£M 


then  the  relative  error  is 


/(^)//nw=  n 

k~N+l 


l  +  (z/ik? 
1  + 


n 


1  +  iz/ijk  j 


X 


n 


l-'izhh) 


Actually,  as  we  shall  see  in  a  moment,  this  is  not  the  best  possible  approximation;  to 
improve  the  convergence  at  infinity  one  needs  one  minor  correction. 

Our  estimates  on  the  error  are  based  on  the  following  result: 
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Lemma  4.1.  For  each  positive  integer  N  define 


Proof,  For  a  proof  for  the  first  part,  and  for  additional  terms  in  the  expansion,  sec 
[14,  p.  34,  formula  (14)],  and  to  get  the  second  part,  use  the  first  part  together  witli 
(28).  D 

Now  let  us  apply  this  to  the  rod  of  constant  unit  density.  In  this  case  we  get 

f,  ...  ,  ,  T{N  +  l/2)^r(iz/ff  +  iV  +  l)r(-iz/7r  +  jV  +  1) 

^  l)2r0^/5r  +  N+  l/2)r(-iz/7r  +  N+1  /2)' 

If  we  expand  the  functions  above  as  ^  co,  2  >  0,  by  using  the  asymptotic  for¬ 
mulae  (29),  we  observe  that  for  large  values  of  2:  the  ratio  beliavcs  like 

(1  +  z^I{'k{N  +  1/4)^)*^^.  Thus,  we  get  a  much  better  approximation  if,  instead  of  us¬ 
ing  the  approximation  that  we  originally  proposed  to  use,  we  divide  by  this  square  root. 
Thus,  the  approximation  that  we  end  up  with  for  C  after  this  first  stage  is 


f{s) 


(31) 


f{s) _ 


S(s)  (1  +  py{7c{N  +  l/4)yf^  Ui  1  -  1/2))' 

which  leads  to  a  relative  error 

1 

^iV  /'I  ,  07//^/ 

(32) 


i+py{irky 


tW(s)  =  =  — _ 

C(s)  (i+^2/(,r(iV  +  l/4))2)'/'“ 


r{N  +  \l2fT{\0lTt  +  .V  +  l)r(-i/?/x  +  N  +  1) 
r{N  +  lfT{[p/%  +  N  +  ll2)r(-iplz+N  +  1/2)' 

for  the  rod  with  density  p  =  1.  The  asymptotic  estimates  in  Lemma  4.2  can  now  Ix' 
used  to  estimate  the  relative  error  and  we  obtain 

Theorem  4.3.  Let  A(t)  be  as  in  Section  3  and  let  be  the  relative  eiror  for  th.e 
rod  with  constant  density  /?(.r)  =  1  defined  in  equation  (32)  above. 

1.  For  s  in  a  compact  subset  K  of  Us  >  0, 


(33) 


4'’(s)  =  1  +  0(Af-2)  asN^oo 


2.  If  p  maps  >  0  into  a  region  that  is  asymptotically  contained  in  a  sector 

|arg2:l  <  I  -  €,  e  >  0,  (recall  that  this  implies  yl(Od-)  =  00),  then 

(34)  =  0(JV-2)  asN-^00. 

3.  If  ->1'(0+)  =  00,  then 

(35)  “  l||//«  “4  0  as  fV  00. 

If,  in  addition,  A(O-l-)  <  00,  then  we  have  the  more  precise  estimate 


(36) 


\\,[«.  =  0  as  N  00. 
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If  —A'{0)  <  oo,  then 

(37)  limsupl|T^^^||/roo  <  oo  and  limsupyi/ 

N—*oo  N-*oo 

In  fact,  ifW{s)  is  any  weight  function  with  |W^(5)|  — >  0  as  |5|  — *  oo  in  3Ra-  >  0. 
then 

(38)  -  0  and  -  1)|1//-  0  as  N  - 

In  addition, 

limsup|4*^(iw)|  (l  +  0{N--))  , 

(39) 

"TTc 

as  N  oo  where 


(40)  C  =  exp  (a'(0+) M(0+)^/^)  . 


Proof  We  use  Lemma  4.2  to  estimate  the  right  hand  side  of  (32)  as  follows. 
First,  formula  (29)  with  c  =  1/2  and  d  =  1  yields 

(41)  ^ry+^i'f'  = 


Since  the  region  11  is  symmetric  with  respect  to  the  real  axis,  it  suffices  to  consid<'r 
the  case  where  s  is  in  the  first  quadrant  Q={5|9f?s>0,Qf5>0}. 

If  s  G  Q,  then  p{s)  G  Q  and  hence  3^^(~i/?/7r+  +  1)  >  N  for  s  G  Q.  Tliu.s  Ijy 

(29)  with  c  =  AT  4- 1  and  d  =  N  1/2, 


r(-i/?A+jv  +  i) 

^  ’  r(-Wir+Ar+l/2) 


{-iPlir  +  N+ll^f^^  X  [l  +  0(A^-^)]  as  A' 


CO, 


where  the  O  term  holds  uniformly  for  s  G  Q. 

It  remains  to  estimate  the  difficult  factor  in  (32),  namely  the  ratio 


r(i/?/7r  +  iV  +  l) 

r(i/?/^  +  iv  +  i/2)' 


If  s  lies  in  a  compact  subset  K  of  Q,  or  if  the  region  FI  is  asymptotically  contained 
in  a  se^jtor  |arg2|  <  |  —  6  for  some  e  >  0,  we  can  again  use  (29)  with  c  =  N  +  1  and 
d  =  AT  +  1/2  to  obtain 


rjW/i^+N  +  i) 

^  ’  T(iPlv+N +  1/2) 


{ip/iT  +  N  +  1/4)'^^  X  [l  +  0(A^-^)]  as  N  oc, 
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for  all  s  6  K,  or  for  all  Q  in  case  the  region  11  is  as}'mptotically  contained  in  |arg  .;|  < 

.2  -  In  the  latter  case  the  constant  in  the  O  terra  depends  on  how  close  e  is  to  zero. 
Combining  (41),  (42)  and  (43),  completes  the  proof  of  Parts  1  and  2  of  Theorem  4.’ , 
When  A'{0+)  =  -oo,  but  the  region  H  is  not  contained  in  any  proper  subsector  ol 
the  right  half-plane,  we  again  use  (29)  with  c=  N +  l,d=  N +  1/2  and  £  -  5r/2  to  see 
that  as  N  oo 

f44i  I!Ml±E±Jl.==(i^/T!  +  N  +  l/if^><  [l-l-0((i/3/irH-A'-l-l/4)--)] 
r(i/?/7r  + AT  4-1/2)  ^  ' 

uniformly  for  those  s  for  which  Rms)/^+N+ 1)  >  0.  If,  on  the  other  hand  3?(i;l/-4- 
At  +  1)  <  0,  we  use  (30)  with  c  =  4- 1,  d  =  A1  4- 1/2  and  £  -  7r/2  to  obtain 

m/^  +  M+1)__  sin 

TiKi/n  +  N  +  1/2)  ^  siu(!/?4-(A  4-l)tr) 

X  [l  4-  0  ((-i/J/tr  -  Af  -  1/4)"")] 

=  (-i/3/ir  -  A'  -  1/4)'^^ e’'/^  (tT^) 

X  [14-0  ((-i^/tr-Af- 1/4)--)] 

=  (i/?/tr-bA4-l/4)‘/' 

X  [l4-0((i;8/tr4-A4-l/4)--)]  , 

as  V  -»  00.  Here  the  branch  cut  has  been  chosen  to  make  (i/?/7r4-A'‘ 4-1/4)'/^  coiitimions 
for  Ks  >  0.  Combining  (44)  and  (45)  with  (41)  and  (42)  we  find  that  for  .s  6  0, 

(46)  =  1  +  0{N-^)  4-  O  ((i/?/ir  4-  4- 1/4)"^)  as  A  oo 

when  »i(i/?/7r4-  A  4- 1)  >  0,  and 

(47)  ^(t)(,)  3=  (i  +  0{N-^-)  +  0  ((i/?/tr  +  N  +  1/4)-^))  as  N  -  oo 

when  S(i/J/tr  4-  A  4- 1)  <  0.  By  Part  3  of  Lemma  3.2,  3l/4(M  -»  oo  as  adiia--)  -  o: 

when -A'(04-)  =  00,  and  it  follows  frob  (46)  and  (47)  thatllr.v  -1||„»  ->0asA  ^  o: 

as  asserted  in  (35).  .  *  .  •  , 

For  the  more  precise  estimate  (36),  note  first  that  Part  1  above  permits  us  to  irstnc 

attention  to  a  subset  Q'  =  Q\  K,  K  compact,  and  use  the  asymptotic  estimate  and 
inequalities  (23)  and  (24)  of  Part  2  of  Lemma  3.2.  By  Part  2,  we  may  also  assunie  tl.ai 
ir/4  <  arg/?(s)  <  tr/2  for  s  £  Q'.  Thus,  for  each  s  €  Q',  there  exists  ai{,s)  €  R  with 

(48)  {iu{s))  =  9iP{'iu>{s)) 

Let  UJK  denote  any  positive  number  with  a/J(ia;,v)  =  By  (23)  and  (24)  il  is 
sufficient  to  prove  (36)  with  K/?(i  A)  on  the  right  side  replaced  by  3J,4( iu,-.v )■ 


For  those  s  in  Q'  with  >  {N  +  l)7r  (i.e.,  where  we  use  (47)),  wc  get  from  (4S) 

and  (23),  (24),  the  uniform  estimate 


It  follows  that 


(50) 


1+e-W 

1  _  e-2«^) 


1  +  O 


as  N  -»  00,  uniformly  for  s  &  Q'  satisfying  Ji(i;9{s)/5r  +  +  1)  <  0. 

Next,  note  that 


(51)  |i/3(s)/tr  +  N+  1/4|  >  max{3}/?(s)/5r,  |-a/?(s)/7r  +  N  +  1/4|} . 
Also  observe  that  by  Part  2  of  Lemma  3.2, 

(52)  3fJ^(iu/>)  =  o{N)  (jV  —V  oo). 

Thus,  unless 


(53)  \-^P{s)/7r  +  AT  +  1/4|  >  3^,i^(iwAr), 
we  will  have  1/2  <  u{s)/ojn  <  2  by  (23),  and  hence,  by  (24), 

(54)  ^  ^ 

Combining  (51),  (53)  and  (54)  gives  us 


(55)  |iy3(s)/5r  +  A^+l/4|-i  =  0{(»/?(iw,v))"*)  (A' oc) 

uniformly  in  Q',  and  (36)  follows  from  (46),  (47),  (50),  (52)  and  (55).  This  conipk’tc.s 
the  proof  of  Part  3  of  Theorem  4.3. 

Finally,  if-A'(0+)  <  oo,  then  Part  3  of  Lemma3.2  gives  Q:/5(iLj)  -  u;A(0h-)-^A'  and 
9f?/?(ia;)  -|A'(0+)/A(0+)^^^  >  0  as  a»  oo.  In  particular,  the  curve  \3{i<jj)/-k  (u;  > 

0)  lies  in  the  second  quadrant  but  remains  outside  some  strip  { 9^::  <  -m,  0  <  < 

-“A'(0+)/4A(0+)^/^};  by  the  argument  principle,  {  i/?{s)/7r  |  s  €  Q}  is  disjoint  from 
this  strip  as  well.  It  follows  that  when  >  m  +  1, 

|i^(s)/7r  + A^  + 1/4|  >  M  (s  €  Q) 

and  that 


1  +  e-2/3(.q 

1  _  e-ipu) 


<M 


foi  those  s  €  Q  with  9?(i^(s)/7rH-A^4'l)  <  0.  By  (46)  and  (47),  we  get  the  first  inequality 
ill  (37).  For  we  apply  the  asymptotics  of  Lemma  4.2  to  the  reciprocal  of  (32): 
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this  leads  to  analogues  of  (46)  and  (47)  from  which  the  bound  on  1/rV  in  (37)  can  ho 
deduced  as  above.  Comihning  (37)  with  (33)  we  see  that  (38)  holds  whenever  IK  is  a 
weight  with  lVV'(s)|  ->  0  as  s  -t  oo  in  5R.s  >  0.  The  maximum  oscillation  estimate  (39) 
is  an  immediate  consequence  of  (47)  and  the  fact  that  R/7(iiv)  — *  — (0+).4(0+) 
as  w  oo.  This  completes  the  proof  of  Part  4  of  Theorem  4.3.  □ 

For  the  particular  viscoelastic  models  Ai  ~  At  in  Section  3,  we  have 
Corollary  4.4.  Let  be  the  relative  error  for  the  rod  with  constant  density 
p[x)  =  1  defined  in  equation  (32).  Then 

L  For  the  Kelvin-Voigt  model  A  ==  Ai  and  for  the  fractional  dcmudxue  viodel 
A  =  A2, 

(56) 

Here  the  constant  in  the  O  estimate  depends  on  whether  one  has  Ai  or  Ao  and. 
in  the  latter  case,  on  the  order  I  —  fi  of  the  fractional  derivative. 

2.  For  the  intermediate  model  A  =  A3  of  order  1  —  /./  , 

(57)  -  III//0C  =  0(./V-20-^‘))  a..  .V  00. 

3.  For  the  standard  linear  solid  A  =  A^, 

limsup|'rj,J\ia;)l  =  (l  +  0(n~^))  , 

(58)  \lc, 

liminf|4)’(iw)|  =t^  (l  +  ^)) 

|u>|-*oo  i  +  O  ' 

as  N  — >  00  where  C  =  exp  4-  . 

Proof.  Pai'ts  (1)  and  (3)  are  immediate  consequences  of  Theorem  4.3  and  th<- 
discussion  following  Lemma  3.3.  The  size  of  the  constant  in  the  O  term  in  (56)  (k'lieiids 
on  the  sector  opening,  that  is,  on  the  order  1  -//  of  the  fractional  derivatir'e.  For  Pait  2. 
we  use  the  fact  that  for  A  =  .43,  ~  C'(3,3(iui))'-''  ~  C(ui)'-"  as  is  noted  in  the 

discussion  following  Lemma  3.3.  □  ^ 

Remark  4.5.  For  the  general  case  of  cojistant  density  fj{x)  =  />,  the  transfer 

function  for  our  rod  problem  is 

p(,)  =  ik^Mcoth(,/p/J(s)). 

In  this  case  the  appropriate  first-sUige  relative  error  r,!,*’  is  defined  as  vn  cquahon 
above,  except  that  each  occurrence  of  P  is  replaced  by  {^,6).  >ln  cxarninaHon  of  the 
proofs  (the  key  difference  is  that  the  term  (l+e-2^)(l-e-2^)-'  infonmda  (47)  becomes 
(1  +  c-^Ke'^)(l  -  e"^'^'’)'’  ),  shows  that  Theorem  4.3  and  Corollaiy  4.4  renvan  true  for 
the  rod  with  p(x)  =  p  with  the  constants  C  now  given  by 

(59)  C  =  exp  (v^2l'(0+)M(0+)’/^)  andC  =  exp[-,/pS((E  +  6er''/'^, 

23 


respectively, 

In  order  to  treat  the  rod  with  general  densities  (as  well  as  the  beam,  sec  [10])  we 
need  to  take  account  of  the  fact  that  the  zeros  and  poles  of  the  transfer  function  arc  not 
spaced  at  exactly  equal  intervals.  Recall  that  for  the  rod  with  general  density  p{x),  the 
transfer  function  P  is  given  by  P{s)  =  /?(s)s-2/(/3(s))  where  /,  defined  in  Proposition 
3.1,  has  zeros  iigt  and  nonzero  poles  iitj*  satisfying 


(60) 


4  = 


(2k  -  l)7r 

2pi 


+  W, 


JCTT 

'nk  =  —-\-  Vk 
Pi 


{k  oo). 


Here  pi  =  p{xy^^dx,  and  the  perturbation  terms  pk  and  satisfy  the  estimate 


(^1)  Pk  =  0{l/kl  Uk  =  0{l/k)  {k  oo). 


We  assume  that  the  f*.  and  7]k  are  known  to  any  degree  of  accuracy  for  I  <  k  <  N,  ami 
wc  define  the  first-stage  approximation  to  P  by 


Pl}\s}  =  ^  (l  +  {Pi0f/{n{N  +  1/4))=)'- -  ii 


i  +  (/?/6)^ 
i  +  (0hhY' 


where^^by  Proposition  3.1,  po  —  /a  p(a;)da;.  Then  the  first-stage  approximate  compeu- 
sator  for  the  optimal  compensator  C(s)  =  f{s)l{S(s)P{s))  is  given  by  Cy'(s)  = 
'P(^)I(P{^)Pn^(^)),  and  the  relative  enor  is 


(62) 


r'^\s)  = - \ _  ft  ^  + 

(1  +  (pi/3)V(^(tV  +  1/4))2)'/2  1  + 


In  order  to  estimate  the  relative  error  T^\s)y  we  set 

(S3)  ft  =  --g-  and  ^  ,  (k  oo), 

^Pl  Pi 

and  write 


(64) 


_ _  fr 

(1  +  iPi0fmN  +  l/4))2)*/^  1  +  (P/iJ^ 

X  ft  l±i^/^x  ft  ^  + 

k-N+l  k=N+l  1  +  (/9/%)-^ 


1  +  iPH)- 


=  f(l)c,’l  TT  '  „  TT  _ 

^  *=Ar+l  1  +  i^/CkY  k=N+i  3  +  (/?/%)^ 


Here  is  given  by  the  expression  on  the  right  hand  side  of  (32),  but  with  eacli 
occurrence  of  /?  replaced  by  (piP).  Theorem  4.3  and  Remark  4.5  apply  to  4'*.  Tlic 
following  Lemma  yields  estimates  for  the  remaining  two  products  in  equation  (64). 
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Lemma  4.6.  Let  ijik,  4>k  {k  =  1,2 — )  satisfy 

■{05)  V-I.  >0,  Wk-ck\<b  <co 

for  fixed  positive  constants  b  and  c.  Write  6k  =  IV’L  ~  V’*l  tmd  assume  that 
(66)  4  =  0(l/k)  {k  ^  00). 


Define  by 

(67) 


0jv{2) 


Then  as  N  ->  oo  the  followiny  uniform  estimates  hold: 

j  |Q^,(r)  -  1|  =  0(Ejt>A^4A^)  as  N  ^  oo  uniformly  in  compact  subsets  A 
of  C 

a.  \Qn{z)  -  1|  =  0(Zk>N  h/k)  asN^oo  uniformly  in  I>,  =  { .r  |  |arg.2l  <  0  < 

TT  } . 

3,  IQ,, (2)  -11=0  (Zk>N  6k/k  +  sxxpk>N  logA-)  N  ^  oc  uniformly  in  D-,  = 
{  z  I  >  0 }  U  {  2  =  .r  +  ij/ 1  a- <  0,  |t/|  >  d  >  0  } . 

4.  \Qi,{z)  -  1|  =  0{supfc>,v6i.fc'’)  as  N  oo  uniformly  in  D.i  =  { 2  I  iR-:  > 
0  }  U  {  z:  =  .r  +  ij/ 1  .i'  <  0,  jj/l  >  d/{\x\’'  +  1)  }  with  0  <  p  <  1. 

We  reraark  that,  for  each  part  of  Lemma  4.6,  the  conclusion  is  valid  provided  only 
that  the  {4}  are  such  that  the  0  term  in  the  conclusion  is  finite;  this  is  clearly  the  ea.si' 
when  (66)  holds.  Part  4  of  Lemma  4.6  is  not  used  for  our  analysis  of  vi.scoelastic  rods 
since  the  region  H  is  asymptotic  to  a  closed  half-plane  strictly  contained  in  >  0.  M  <• 
include  it  here  since  its  proof  is  a  minor  modification  of  that  of  Part  3,  and  we  will  u.se 
Part  4  to  study  the  analogous  problem  for  a  viscoelastic  Euler-Bcrnoulli  beam  [10]  wth 
6k  -  0(e-'')  as  A-  ->  oo.  This  is  needed  due  to  the  fact  that  for  the  case  of  a  standard 
lineal-  solid  Euler-Bernoulli  beam,  the  image  of  Ss  >  0  under  the  appropriate  analogue 
of  (3  approaches  the  axis  at  infinity  like  a  hyperbola. 

Proof  Throughout  the  proof  M  denotes  a  positive  constant  who.se  value  may 

change  from  one  line  to  the  next. 

Write 

logQ^iz)=  £  iogil +  tk{z)), 

k=N+l 


where 


Since 


*11-  ^  ^  i>'k  -'kk  2 

*  ”  1  -1-  z/xl)[  xpk  z  d-  ip'k 


|Qiv-l|  < 


/: 


logQw 


<  llogQ,v|exp|logQA'l 
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(the  integral  is  the  path  integral  along  the  straight  line  joining  0  to  log^Ar),  it  is  clearly 
sufficient  to  show  that  the  conclusions  of  Parts  2,  3  and  3  hold  with  their  loft  sides 
replaced  by 


Tn{z)  = 

h=N+l 

But  in  Z)i,  \z{z  +  is  uniformly  bounded,  so 

(68)  T4z)<M  f;  ^  (3  6  A). 

k=N+l 

since  ipk  ~  cfc,  and  Part  2  is  proved.  For  2  in  a  compact  subset  K  of  C,  [3(3  +  '  |  = 

0(k~^)  uniformly  as  -»  00,  and  we  obtain  an  extra  factor  of  k  in  the  denominator  of 
the  sum  in  (68).  This  proves  Part  1. 

To  prove  Part  3,  we  can  restrict  attention  to  D'o  =  £»2\£»i  (9  =  Sw/i).  For  2  =  x+iy 

in  D'2,  if  if>i  <  -x/2  or  >  -3a/2,  then  |z(z  +  <  M,  and  we  can  proceed  as 

above.  Let 

A  w(2)  -{k  I  k>  N  and  -  x/2  <f'^<  -3x/2  },  Ak(z)  =  sup  i\. 

A;€A’,v(=) 

(with  A^r{z)  =  0  if  Ki^{z)  is  empty).  Then 

I)v{^)=  E  Uz)\<  M  E 

(69)  A-€A',v(c)  ka<N{:) 

<MAs(z)  y  r-i — 

Ma<N{z) 

since  \zli)k\  is  uniformly  bounded  above  and  below  on  A'.v(^). 

To  estimate  the  sum  in  (69),  wc  use 

\^  +  ‘<P'k\  ^  niax{|y|.  + 

>  max{|2/|,  |a-  +  cfc|-6}. 

For  fixed  z,  if  Aa-(2)  >  0,  there  is  aj  e  Kfi{z)  such  that  6j  >  iA.v(2).  In  particular. 
j  >  N  and 


cj  +  b>  —  >  — ^Izl. 

2  2V2 


It  follows  that 


E  I  ■  r/. 

-tCAVW  l~+^l 

f 

<2Sj  y  y  - 1 _ 

I  teAwW  1^1  teAX--)  “  *) 

\(x+cfc|<6+c  |a:-fcfc[>6+c 


The  index  set  for  the  first  sum  contains  no  more  than  2(</c  +  3  terms,  while  the  .s('coiid 
sum  is  comparable  to  the  sum  of  two  finite  harmonic  series.  We  get 

(70)  '  ' 

<M6dl+\ogj)  <  M sup{41og/r}; 

this  proves  Part  3.  The  proof  of  Part  4  is  identical,  except  that  (70)  becomes 

T;.(r)  <  M6j(lt:|'’  +  l+log|j|) 

<  M  sup{6itA:''}. 

k>N 

This  completes  the  proof  of  Lemma  4.6.  □ 

Clearly,  Part  4  can  be  generalized  from  D3  to  other  regions  where  \y\  >  d/(l  +  ) 

is  replaced  by  \!j\  >  ii{\x\)  with  ^i{\x\)  -*  0  as  [.I'l  -»  00. 

Combining  Lemma  4.6  with  Theorem  4.3  and  Remark  4.5  tve  obtain 
Theorem  4.7.  Let  A(t)  and  p{x)  be  as  in  Section  S  and  let  ry  be  the  first- stayc 
relative  error  defined  in  equation  (62)  for  the  rod  with  density  p{x).  Then  Parts  l-f  to 
Theorem  f.S  remain  valid  provided  that  the  following  modifications  are  made. 

1.  In  Part  2  the  estimate  (3f)  must  be  replaced  by  the  weaker  estmuiU 

(71)  -  l|l;/»  =  0(N-')  as  N  00. 

2.  In  Part  3  the  precise  estimate  (36)  tlud  holds  when  .4(0+)  <  00  must  he  replacal 
by 

(72)  114"  -  1|1;,~  =  O  ((S},6(i.V))"")  +  0  (A'-'  logA')  as  oc.  , 

3.  In  Part  4  equation  (39)  must  be  replaced  by 

limsup|rv'(iu))l  ^  ^  (l  +  0(.hi  '  log  .A'))  1 

US)  \-c, 

lira  ^  (l  +  0[N  ^  log  AO) 

|w|— 'OO  i  “T  t-' 

where  the  constant  C  is  now  given  by 

(74)  C  =  exp(piW(0+)/.4(0+)’/4 


wit/ipi  = /o4W’^^d.r.  .  .  .  , 

Proof.  As  noted  earlier,  satisfies  the  conclusions  of  Tlieorcui  4.3  with  loniiul;\ 
(40)  replaced  by  (74).  The  remaining  two  products  in  (64)  may  bo  rcwritt.iMi  as  lour 
products  with  linear  Mobius  factors  each  of  which  may  he  estimated  hy  Lonima  4,G. 
Specifically,  one  of  these  products  is 


Zfii{s) 


■TT  ^ 


Since  (61)  holds  and,  Lemma  3.2,  Hi  z=  {  i/?(s)  |  >  0  }  is  always  nsjanptorically 

contained  in  a  region  of  the  form  D2,  it  follows  from  Part  3  of  Lemma  4.6  that 

l|Zv  ~  l||//oo  =0  f  4.  snp^ifclogA: 

\k>N  ^ 

-0  I  53  +  sup  ^-logA: )  =0  (N~'^  log  A')  as  (A'  oc). 

\k>N  ^‘>^V  J  ^  ^  ^ 

If  n  is  asymptotically  contained  in  a  sector  {z  |  |z|  <  i9  <  7r/2  },  then  Part  2  of  Lemma 
4.6  yields  \\Z^r  -  l||//«>  r=  0{A/’~'^),  while  Part  1  shows  that  Zi^r{s)  =  1  +  0{N~-)  as 
N  00  uniformly  for  s  in  a  compact  subsets  K  of  >  0.  The  other  three  products  of 
linear  Mobius  factors  which  occur  in  (64)  may  be  treated  in  exactly  the  same  maiiuor, 
and  Theorem  4.7  follows.  □ 

Of  course,  the  general  density  version  of  Corollary  4.4  for  the  particular  viscoela.stic 
models  Ai  —  also  holds  after  making  the  same  modifications  as  in  Tlieorem  4.7. 

For  each  of  the  models,  Figure  2  shows  |rj^?\ia;)  -  1|  (ordinate),  as  defined  in  (32). 
plotted  against  u  for  N  increasing  from  N  -  o  to  N  =  55  in  increment.s  of  10.  The 
curves  move  down  and  to  the  right  as  A^  increases. 
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5.  Approximation  of  Mobius  Transforms  in  /?*  by  Mobius  Ti'ansforms  in 
s.  After  the  first  approximation  step,  we  arrive  at  a  finite  product  of  factois  (Mobius 
transformations  in  /?^)  of  the  type 


i+(/g/a) 
i  +  Whk? 


Substituting  the  definition  of  0  into  these  factors  we  get 


Mk(s)  = 


s/ek  +  A{B) 

si4+A(sy 


How  do  we  approximate  these  factors  rational  functions  of  6‘? 

Let  us  take  a  closer  look  at  one  of  these  factors,  where  for  simplicity  we  donotc 
by  f  and  by  ?/,  i.e.,  we  write  it  in  the  form 


(75) 


M{s) 


s/e Ms) 

s/rp-^-AisY 


As  s  00  this  factor  tends  to  which  is  close  to  1  (since  i/He  tei^ds  to  1  as 

k  oo),  and  as  s  0  it  tends  to  1  (since  A(s)  ~  E/s  as  0).  The  largest 
deviation  from  Ion  the  imaginary  axis  occurs  in  the  cross-over  region  where  the  order 
of  magnitude  of  A  is  the  same  as  the  order  of  magnitude  of  s/e  and  s/if.  This  indicates 
that  a  rational  approximation  of  M(s)  should  have  its  zeros  and  polos  somewlu^re  in 
the  left  half-plane  close  to  the  cross-over  region.  As  we  shall  see  below,  this  is  true,  at 
least  if  we  restrict  our  attention  to  the  case  where  the  internal  damping  is  weak  enough. 

In  Section  4  the  amount  of  internal  damping,  described  by  the  parameter  f  in 
our  examples  A1-A4,  was  significant  only  in  that  it  affects  the  size  of  the  O  constants 
in  our  estimates,  and  hence  the  number  of  terms  N  needed  to  get  a  good  fir.st-stage 
approximation.  It  did  not  determine  the  shape  of  the  image  of  the  half-plane  >  0 
under  the  map  /?,  and  consequently  did  not  affect  the  order  of  our  convergence  estimates. 
In  this  Section  we  shall  make  use  of  the  fact  that  in  mam'  materials  of  interest,  tin' 
damping  parameter  e  is  quite  small  compared  to  the  elastic  parameter  E. 

More  preciselj',  we  assume  that  A  is  of  the  form 


(76) 


A(s)  =  E/s-^€d(s), 


wdiere  a  has  been  scaled  so  that,  for  example, 


a(0)  =  1. 

Then  (75)  can  be  written  in  the  form 


(77) 


,  ,  ^  s/^^  +  E/s  +  €a{s) 

^  ’  s/if  +  E/s +  (a(sy 


For  6  =  0  (no  damping)  this  function  has  two  purely  imaginary  zeros  located  at  ±i  \/E^. 
and  two  purely  imaginary  poles  located  at  ±i  y/Ei].  By  continuity,  for  small  uonzcTo 
values  of  6,  M  Avill  have  zeros  and  poles  close  to  these. 
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As  the  following  lemma  shows,  these  are  the  only  complex  zeros  and  polos  of  the 

function  M:  .  o  •?/  x  i  . 

Lemma  5.1.  For  each  constant  C  >  0,  the  function  s/C'  +  A{«) 

pair'  of  complex  conjugate  zeros. 

For  a  proof,  see  Desch  and  Grimmer  [6]. 

For  each  C  and  €  such  that  s/C^+£:/s+e3(s)  has  a  pair  of  complex  conjugate  root.s. 
denote  the  root  of  this  function  in  the  upper  half-plane  by  Then  the  appio])iiat(' 
first  order  approximation  of  M{s)  is 


{1  —  s/s,|,<:)(l  ~ 

(We  have  normalized  N  so  that  A'(0)  =  M{0)  =  1.)  Let  N,(s)  bo  the  same  function 
with  ^  replaced  by  and  replaced  by  Then  the  total  relative  error  that  w<. 
introduce  at  this  stage  is 

The  analytic  estimates  that  we  are  able  to  prove  for  the  size  of  rv’  are  still  in¬ 
complete  \\'e  have  been  able  to  show  that,  under  quite  general  assumptions,  the 
F/“-uorra  of  each  factor  M,(s)/Nt{s)  -  1  is  of  order  0(e)  as  e  ^  0.  Hmvever,  the 
0(€)-constant  that  we  are  able  to  obtain  deteriorates  as  one  multiplies  h  sucee.ssiv.' 
factors  and  lets  N  -*  oo.  This  leaves  open  the  question  of  whether  it  is  in  fact  iriic 
that  ||T,!,r>  -  ll|„~  =  0(f)  as  e  0,  uniformly  in  N.  Obseiwe  that  for  tlie  Kelvin-\oigl 
model  Ai  each  factor  M,(s)  can  be  expressed  as  the  ratio  of  quadratics  in  .v.  so  the 
result  is  exact;, i.e.,  there  is  no  error  introduced  at  this  stage.  W'c  are  also  able  to  .sliou 
that  if  «(<)  is  a  finite  sum  of  exponentially  decaying  terms  (the  natural  generalization 
of  the  standard  linear  solid  model  A^)  so  that  A(s)  is  a  rational  function  (in  which  ease 
the  approximation  of  the  factors  M,{s)  by  N,{s)  serves  only  to  lower  the  order  ol  the 
approximate  compensator),  then  it  is  indeed  true  that  ||r,v  —  ll|;;~  —  0(e)  as  e 
uniformly  in  N.  For  the  models  A-2  and  Ai  we  arc  forced  to  simply  compute  the  error 

uumcricall}'.  ,  r  11 

III  Figure  3,  "  M  (ordinate)  is  plotted  against  u;  for  cadi  of  the  iiiodc'ls 

do  Ai  and  4|.  In  each  case,  N  increases  from  5  to  55  in  increments  of  10  from  tin- 
bo'ttom  graph  to  the  top  one.  Similar  plots,  not  shown  here,  were  made  for  a  iiiore 
refined  procedure  where  a  real  third  root  of  the  numerator  and  of  the  denominator 
of  each  Mk(s)  (these  approach  the  negative  real  zero  of  .4(s)  closest  to  the  origin  as 
k  -»  oo)  was  included  in  the  approximation  of  Mk(s)  to  account  partially  loi  the  tieep 
response.  For  the  parameters  used  here  improvement  by  up  to  a  factor  of  2  iii  the 
relative  error  was  obscn'cd  in  the  irrational  cases  (and  there  is  then  no  error  at  this 
stage  for  model  /I4). 
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6.  The  Final  Approximation.  When  we  combine  the  different  expansions  de¬ 
scribed  in  the  previous  sections  we  get  a  rational  approximation,  except  for  tlie  fact 
that  there  is  a  left-over  irrational  factor  in  the  denominator  of  (31)  that  has  not  Ix-cii 
accounted  for,  namely  the  factor 

(78)  + 

Here  we  look  at  low  order  approximations  of  (78).  In  particular,  we  investigate  the  siz(' 
of  the  error  that  one  gets  by  ignoring  this  factor  completely;  note  that  it  tends  to  1 
uniformly  on  compact  subsets  of  the  right  half-plane  as  N  —>  oo. 

The  relative  errors  rl?’  (jV  =1,2,...),  together  with  their  reciprocals,  arc  umlormly 
lioundcd  in  (Theorems  4.3  and  4.7),  and  the  same  is  true  of  the  errors  r^,  of  Sect  ion 
5,  at  least  in  the  special  case  whore  A(s)  =  E/s+ea(s)  is  rational  and  e  is  small  rolahtc 
to  B  This  type  of  boundedness  will  no  longer  hold  for  the  new  error  introduciid  when 
the  square  root  in  (78)  is  approximated  by  a  rational  function,  because  in  general  Hus 
root  will  not  have  a  rational  growth  rate  at  infinity.  Thus,  in  the  final  approxunatiou. 
in  most  cases  the  magnitude  of  must  be  either  unbounded,  or  not  boundetl  away 
from  zero.  However,  as  we  saw  in  Section  2,  the  only  thing  that  wc  really  have  to  worry 
about  is  that  we  satisfy  the  requirement  that  limsup,,’_coll'^-vC51|;/»^  <  oc-  More 
precisely  since  has  been  unifonnly  bounded  in  the  preceding  sections,  if  wo  donot.- 
thc  approximation  that  we  will  use  for  the  root  in  (78)  by  /i,v,  then  we  need  only  harc 
h.v  converging  to  one  uniformly  on  compact  subsets  of  the  right,  halt-plane,  and.  m 
addition, 

(79) 


where 


CAs)  a  +  l3{s)-/MN  +  l/^)?)''-  ■ 
C'-^\s)  ~  'bv(s) 


At  this  stage  the  size  of  the  function  C  at  infinity  becomes  important.  There  are  iw<> 
pos.sibilitiGs:  either  C(oo)  =  0  or  limsup,v_ool<5(s)l  >  0.  (In  the  former  case  the  ideal 
compensator  is  .strictly  proper,  in  the.  latter  it  is  not.) 

The  second  of  these  possibilities  can  occur  in  only  one  way.  Recall  that  c  (.s)  ~ 
f(s)/P(s)S{s)  --  f(s)/P{s)  as  s  00.  The  function  f  must  have  a  zero  of  inU'j-vi- 
order  at  infinity,  since  we  require  it  to  be  rational.  Concerning  P{s).  wc  have  th(‘ 
following  estimates  from  the  results  of  Section  3. 

Lemma  6.1.  LetA(t)  be  as  in  Section  3. 

1.  P{s)  is  analytic  in  the  closed  right  half-plane,  except  for  a  double  pole  at  0.  and 

1  +  I/5WI 


is  bounded  hnd  botinded  away  from  0  in  >0.  Moreover, 
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2.  if  y4{0+)  <  cx),  then  |6'P{5)|  is  bounded  and  bounded  away  from  0  as  s  ^  oo  in 
9f?5  >  0,  and 

3.  ifA{0-\-)  =  00,  then  |P{s)|  =  o{|5|”^)  and  |P~^(5)|  =  O  as  s  ^  oc  in 

>0.  ^  ^ 

Proof.  By  Lemma  4.6,  we  need  only  look  at  (17),  the  case.of  uniform  dcu.sity.  Then 
Fait  1  follows  from  Lemma  3.2,  Parts  (lb)  and  3  (with  the  argument  principle).  Parts 
2  and  3  are  cleai^from  Lemma  3.2,  Parts  (Ig)and  (Ih).  D 

Thus,  since  C{s)  must  be  bounded  for  large  s  (see  (4)),  the  only  case  where  C'(oo)  = 
0  fails  to  hold  is  the  one  wdiere  -4(0+)  <  oo  and  T  has  a  first  order  zero  at  infinity.  In 
all  other  cases  f  must  have  at  least  a  second  order  zero  at  infinity,  and  C(oo)  =  0. 

Let  us  first  discuss  the  case  that  in  view  of  the  discussion  above  is  more  common, 
namely  the  one  where  C(oo)  =  0,  and  f  has  at  least  a  second  order  zero  at  infinity' 
^Ve  claim  that  in  this  case  \ve  may  take  /ijv  =  1,  i.e.,  w-e  may  ignore  the  square  root 
c.ompletcl}^  This  follows  from  the  follo'wing  lemma: 

Lemma  6.2.  Iff(s)  =  0{$~-)  as  s  ^  oo,  then  (70)  holds  with  =  1. 

Pioof  By  Part  1  of  Lemma  6.1,  it  suffices  to  sho^v  that 


lim  sup  sup 

A’-^oo  3fj5>0 


|l+;97(ff(AUl/4))2|'/2 

1  +  1)9(71 


<  00. 


But  this  is  triviall}^  true  since  the  numerator  can  be  estimated  from  abo\’e  In-  1  + 

|/?(5)l/(7r(A'  +  l/4))<l  +  |/?(5)|.0 

Note  that  if  \ve  use  =  1,  then  the  approximate  compensators  Cy  = 
will  not  be  strictly  proper  w'hcn  T  has  a  second  order  zero  at  infinity.  If  ouc'  wants  a 
strictly  proper  compensator,  then  one  may  instead  use  /?,y(s)  =  (1  +  e.Y-s),  whore  f  v  0 
as  N  oo;  this  sequence  can  be  chosen  more  or  less  arbitrarilv.  This  does  not  disturb 
(79). 

Noav  let  us  return  to  the  case  where  0(00)  ^  0.  Clearly,  in  this  case  we  cannot 
take  h^T  =  1,  since  the  square  root  in  (78)  grows  like  5  as  s  00,  invalidating  (79). 
However,  the  following  modification  is  sufficient. 

Lemma  6.3.  Iff{$)  has  a  first  order  zero  at  infinity  (and  hence  C{oc)  ^  0 ).  then 
(70)  holds  with  hy  =  1  +  e^s,  where,  e.g.,  =  l/(7rA^  +  1/4),  or  more  generally, 

represents  any  sequence  satisfying  ->  0  and  =  0(N)  as  N  00. 

Proof  In  this  case  C(oo)5(oo)  ^  0,  and  ^ve  must  show'  that 


A— too  |1  +  €A^^j 


This  is  true,  because  |1  +  £/,rs|  >  (1  +  f;v|s|)/\/2,  and,  for  some  constant  A'  >  1. 
and  (ff  >  l/(A'At),  and  hence 


|1  +  pyjitjN  +  1/4))^|‘^^  ^  1  +  A>|/(zAf  +  1/4) 
|1+«H  “  V^(l  +  |.sl/(/i'A'))  ' 

the  supremum  of  which  over  s  stays  bounded  as  A^  — ^  oo.  □ 
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Fig.  4.  Plots  0/ |r5?'(iu;)  -  ![ 


Note  that,  in  contrast  to  the  error  introduced  in  Section  5,  r.v  (.s)  <ii>i>ioa<  lie  s 
one  as  N  -  oo;  there  is  no  residual  error.  Nonetheless,  at  Ica.st  in  our  cxaini-lcs.  the 
error  r},?’  will  be  the  dominating  error  for  high  frequencies. 

We  plot  |tv ’(iu>)  -  1|  (for  N  between  5  and  55,  increasing  to  the  right )  versus  in 

Ficni'G  4.  ... 

To  produce  some  pictures  of  the  final  result  wc  have  fixed  the  ideal  stiisitiMt> 

function  S  by  choosing 


5(6- 


b{s  +  e) 
s2(l  +  es) 


where  6  =  5  and  e  =  0.01.  Then  f  =  1  -  5  has  a  second  order /cto  at  infinity.  Here  .- 
determines  the  cross-over  frequency,  and  e  restricts  the  size  of  PS  at  zero  and  the  sizy 
of  CS  at  infinity;  cf.  (4).  A  plot  of  |51  and  ir|  is  given  in  Figure  o,  as  well  as  plots  ol  C  | 
for  different  choices  of  kernels.  Plots  of  \PS\  and  |CSi,  which  pertain  to  the  sty. b. hi y 
requirement  (4),  axe  given  in  Figure  6.  The  sizes  of  |SCta.1  (see  (13))  plotted  iii 
Figure  7,  and  the  error  |rw  -  Ij  weighted  by  \T\  is  plotted  in  Figure  8.  Here  r.v  i.s  the 
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Fig.  5.  Plots  o/|5|,  lf|  (left,  S(oo)  =  1)  and  ]C\ 


Fig.  C.  Plots  of\CS\  (left)  and  |P5| 

total  relative  error  ta-  =  r»VirV<?\  (As  above,  the  a;  axis  is  horizontal.) 

The_se  plots  summarize  the  combme_d_effect  of  all  the  approximations.  Concernins 
Figure  7,  note  that  the  weighting  by  CS  controls  the  unbounclcclnes.s  of  r.y  at  hi<d” 
frequencies  (due  to  t,v  )  to  the  extent  that  the  validity  of  (13)  is  determined  at  the  fii^l 
peak  of  the  graph  [u  ~  1.5).  In  Figure  8,  obseivc  that  for  Aj,  Aj  and  A.j,  the  curve.s 
rise  to  a  limit  in  the  lower  frequencies  as  N  increases;  this  residual  error  rellec.t.s  the 
neglected  creep  response  in  the  approximation  of  Section  5.  For  A, ,  the  Moliius  factoi-s 
of  Section  5  can  be  expressed  as  quotients  of  second  degree  polynomials  in  .v,  so  no 
approximation  is  necessaiy  arid  we  see  no  “folding  owr”  in  the  figure.  For  .4.1,  wi'  got 
cubic  polynomials;  by  omitting  the  approximation  step  of  Section  5  we  end  up  with  a 
giapli  (not  shown)  similar  to  the  one  for  Ai  above  at  low  frcquencios. 
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Fig.  7.  Plots  of  1t,vC'S|  for  the  final  (rational)  approximation  (o  <  N  <  oo) 
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